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Abstract — A  framework  for  geometric  regularization  of  eleva¬ 
tion  maps  is  introduced  in  this  paper.  The  framework  takes  into 
account  errors  in  the  data,  which  form  part  of  standard  elevation 
maps  specifications,  as  well  as  possible  additional  user/application 
dependent  constraints.  The  algorithm  is  based  on  adapting  the 
theory  of  geometric  active  surfaces  to  the  problem  of  regular¬ 
izing  elevation  maps.  We  present  the  underlying  concepts  and 
numerical  experiments  showing  the  effectiveness  and  potential 
of  this  theory. 

Index  terms:  Elevation  data,  geometric  active  surfaces,  regular¬ 
ization,  absolute  errors,  constraints. 

I.  Introduction 

Elevation  data  often  comes  with  information 
about  measurements  uncertainty.  In  particular,  the 
standard  specification  of  the  National-Geospatial 
Intelligence  Agency  (formerly  NIMA)  for  digital 
terrain  elevation  data  (DTED)  includes  information 
about  absolute  vertical  and  horizontal  errors  (see  for 
example  http://www.fas.org/irp/program/core/dted.htm  and 
http  ://mountains .  ece .  umn .  edu/~  guille/dtedspecification  .pdf 
for  a  copy  of  the  official  DTED  specification.).  Any  given 
post  reading  (e.g.,  with  posts  spaced  approximately  30  meters 
apart  for  the  DTED-2  data  used  for  this  paper)  can  actually 
“move”  up/down  or  to  the  sides,  following  the  absolute 
vertical  and  horizontal  errors  specified  within  the  dataset. 
These  errors/adjustments  can  be,  for  example,  of  5-10  meters 
vertically  and  10-15  meters  horizontally. 

Due  to  such  intrinsic  uncertainties  in  elevation  data,  its 
regularization  is  often  needed,  both  for  visualization  and  for 
processing  (e.g.,  compression).  In  this  paper  we  introduce  a 
framework  to  perform  this  regularization,  while  maintaining 
the  constraints  given  by  the  provided  information  about  verti¬ 
cal  and  horizontal  absolute  errors.  These  errors  are  represented 
in  our  work  by  “uncertainty  cylinders,”  centered  at  the  given 
data  points,  see  Figure  1.  The  cylinder’s  height  represents  the 
absolute  vertical  accuracy,  while  base  diameter  corresponds 
to  horizontal  accuracy.1  We  search  for  a  smooth  surface  that 
deforms  the  original  data  while  keeping  the  points  inside  their 
corresponding  cylinders.  Such  an  important  constraint,  that  is 
dictated  by  the  elevation  dataset  specifications,  has  not  been 
exploited  by  previous  state-of-the-art  approaches  such  as  those 

'if  other  noise  models  are  considered,  e.g.,  Gaussian  distributions,  the  cylin¬ 
ders  are  to  be  replaced  by  other  three  dimensional  shapes  such  as  ellipsoids. 
We  opt  to  work  with  cylinders  since  the  standard  DTED-2  specification  does 
provide  absolute  errors. 


in  [2],  [13].  We  thus  examine  techniques  based  on  geometric 
partial  differential  equations  [8],  [9],  [10],  and  in  particular 
the  geodesic  active  contours  and  minimal  surfaces  framework 
[1],  [5].  In  this  note  we  propose  a  means  to  extend  these 
works  to  the  problem  of  regularization  with  constraints  for 
elevation  maps.  We  show  that  the  framework  is  valid  not 
only  in  order  to  respect  the  absolute  vertical  and  horizontal 
errors  provided  by  the  DTED  specification,  but  can  also  deal 
with  other  constraints,  such  as  the  preservation  of  geometric 
features. 

The  problem  of  regularization  with  constraints  has  been 
addressed  in  the  literature.  Kimeldorf  and  Wahba,  [6],  showed 
how  to  compute  one  dimensional  splines  with  hard  vertical- 
error  constraints.  While  this  elegant  approach  can  easily  be 
extended  to  higher  dimensions,  it  doesn’t  include  the  hori¬ 
zontal  freedom  given  by  the  horizontal  absolute  error.  It  is 
also  not  developed  for  the  additional  geometric  constraints 
that  are  natural  to  add  in  our  framework.  The  theory  of  total 
least  squares  [4]  also  addresses  the  “freedom  of  motion”  of 
the  given  data,  both  in  the  vertical  and  horizontal  position.  In 
it’s  original  form,  although  computationally  very  efficient,  the 
framework  does  not  provide  hard  constraints  (that  is,  the  error 
is  not  guaranteed  to  be  bellow  the  allowed  bounds),  neither 
does  include  any  kind  of  explicit  regularization  or  geometric 
terms.  In  order  to  add  these  important  constraints,  the  problem 
has  to  be  transformed  into  a  variational  framework  much  of 
the  flavor  here  introduced. 

Before  proceeding  with  the  framework  description,  we 
should  comment  that  we  opt  in  this  work  to  regularize  the 
whole  data  set  at  once  (as  a  function  on  the  plane),  and  not 
just  individual  iso-height  lines.  Although  we  could  use  the 
same  framework  for  regularizing  lines,  this  is  a  non-consistent 
approach.  That  is,  when  regularizing  “on  demand”  individual 
iso-height  lines  with  the  absolute  errors  constraints,  terrain 
points  can  be  assigned  to  multiple  elevations.  For  example, 
considering  a  vertical  error  of  10  meters,  a  point  measured  at 
95  meters  could  be  assigned  both  to  lines  at  90  or  100  meters 
if  the  corresponding  iso-level  lines  are  treated  independently. 

II.  Regularization  with  constraints  via 

GEOMETRIC  ACTIVE  SURFACES 

Let  us  begin  with  the  general  description  of  the  minimal 
surfaces  approach  [1],  which  we  then  adapt  to  the  problem  at 
hand.  Let  S  be  a  two  dimensional  (2D)  surface  embedded  in 
the  Euclidean  space  R3.  Let  also  /  :  R3  — >  R  be  some  given 
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Fig.  1.  Problem  description  (data  in  meters):  The  given  surface  is  allowed  to 
move  inside  the  uncertainty  cylinders,  which  are  determined  by  the  absolute 
horizontal  and  vertical  accuracy  of  the  measurements  (as  described  in  the  NGA 
DTED  Performance  Specification  Document  and  provided  in  the  standard 
DTED  file). 


data,  which  constraints  the  surface  (e.g.,  a  gradient  map  if  we 
are  concerned  with  segmentation,  the  original  motivation  for 
the  minimal  surfaces  model),  and  g  :  R  — »  R  a  function  that 
acts  on  this  data.  The  basic  idea  is  to  find  a  surface  S  such 
that  it  minimizes  the  energy  given  by 

J  g(I(S))da,  (1) 

where  da  is  the  classical  element  of  area  on  the  surface  S.  In 
words,  the  above  energy  is  minimizing  the  measure  g  on  the 
surface  S.  In  the  case  of  image  segmentation  for  example,  g 
is  often  a  decreasing  function  of  the  image  gradient,  and  in 
such  a  way,  S  is  constrained  to  be  at  regions  of  large  gradient, 
thereby  being  attracted  to  object  boundaries  in  simple  models. 

The  above  energy  is  often  minimized  via  a  gradient  descent 
flow,  leading  to  the  evolution  equation 

^  =  («5(/($))  -  Vg  •  N)N ,  (2) 

where  k  is  the  surface  mean  curvature  and  N  is  the  surface 
unit  normal.  The  above  equation  is  iteratively  run  to  steady 
state  starting  from  some  initial  guess  So .  It  has  been  shown  in 
the  literature,  [1],  that  the  flow  is  regularizing  the  surface  (due 
to  the  curvature  component)  while  being  constraint  to  move  it 
to  areas  of  low  g. 

With  respect  to  our  surface  elevation  data,  we  select  g  to 
be  a  function  of  the  distance  to  the  uncertainty  cylinders  (see 
also  [14]  for  the  use  of  distance  constraints  in  the  problem  of 
surface  reconstruction  from  scattered  data).  That  is,  inside  the 
cylinders  we  select  g  =  0  (or  some  very  small  e  >  0  value  for 
additional  smoothness),  and  outside  of  the  cylinders  we  select 

g(x,y,z)  =  f(d(x,  y,  z)), 


Another  important  adaptation  to  our  problem  at  hand  is 
that  the  elevation  data  is  not  a  generic  surface  but  a  height 
function  h(x,y)  :  R 2  — ►  R+  (S  =  (x,y,  h(x,y))).  Thereby, 
considering  the  vector  z  a  unit  vector  in  the  vertical  direction, 
the  above  flow  becomes  a  simple  scalar  motion: 

dh  ->  -> 

—  =  (nd(x,  y,  h(x ,  y))  -  Vd(x,  y,  h(x,  y))  •  N) N  ■  z.  (3) 

Here,  the  generic  g  function  has  been  replaced  with  the 
cylinder/distance-based  version  described  above. 

In  addition,  for  functions,  it  is  well  known  that  (subscripts 
indicate  derivatives) 


and 


making  the  numerical  implementation  of  the  above  mentioned 
flow  quite  straightforward. 


A.  Extensions 

The  framework  introduced  above  is  flexible  and  natural  for 
extensions  that  are  relevant  for  the  regularization  of  elevation 
maps.  For  example,  the  following  procedures  may  be  imple¬ 
mented: 

1)  If,  in-spite  of  the  existence  of  the  absolute  vertical 
and  horizontal  errors,  a  certain  feature  in  the  original 
elevation  data  is  not  to  be  allowed  to  move,  all  that 
needs  to  be  done  is  to  place  a  zero  height  and  zero 
width  (or  as  small  as  desired)  cylinder  for  that  point.  If 
for  example  a  ridge  is  to  be  preserved,  after  it  has  been 
detected  by  a  desired  manual  or  automatic  technique,  the 
corresponding  points  in  the  DTED  file  are  marked  and 
their  cylinders  are  reduced  to  the  desired  size  (which  can 
be  null  for  absolute  no  motion).  Note  that  this  can  be 
done  for  a  group  of  points  representing  lines  of  interest 
as  well  as  for  entire  regions.  Of  course,  the  algorithm 
works  for  cylinders  of  different  sizes  and  also  general 
shapes  of  “allowed  motion.”  This  is  easily  handled  by 
the  distance  function  and  from  it  the  penalty  g  function. 

2)  If  other  characteristics  of  the  data  are  to  be  preserved, 
these  can  be  incorporated  in  the  form  of  hard  or  soft 
constraints  to  the  variational  formulation.  For  example, 
a  slope  preserving  term  can  be  easily  added,  leading  to 
an  energy  of  the  form 

J  g(h)da  +  A  J  ||  V/i  -  Vh0  ||2  dxdy, 


where  /(•)  is  any  increasing  function  (in  our  examples  we 
select  the  unit  constant  function),  and  d(x,  y ,  z)  is  the  Eu¬ 
clidean  distance  between  the  3D  point  (x,y,z)  and  the  set  Ci 
of  cylinders  (see  Figure  1). 


where  ho  is  the  original  (noisy)  data  and  A  is  a  standard 
Lagrange  multiplier.  Proximity  constraints  to  desired 
coordinates  can  be  added  as  well  to  force  regions  to 
be  either  preserved  or  moved  towards  a  certain  target. 
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III.  Experimental  Results 

We  experimented  with  DTED-2  data  (30  meters  post  separa¬ 
tion)  with  10  meters  of  vertical  absolute  error  and  13  meters  of 
horizontal  absolute  error,  as  dictated  by  the  data  specification. 
To  capture  the  size  of  these  errors,  we  constructed  a  refined 
grid,  7  times  larger  in  each  horizontal  dimension  (that  is,  the 
original  data  is  linearly  interpolated  to  posts  at  30/7  meters 
to  produce  ho,  the  initial  condition  for  Equation  (3)).  A 
finer  refinement,  which  will  increase  the  computational  cost, 
was  not  found  to  be  needed.  This  data  was  then  regularized 
with  the  technique  described  above,  using  standard  numerical 
techniques  [10].  If  desired,  fast  numerical  implementations  fol¬ 
lowing  [3],  [7]  could  be  used  as  well.  Figures  2  and  3  present 
two  examples  showing  original  and  regularized  results.  Figure 
4  shows  the  histogram  of  the  vertical  motion  (adjustment) 
introduced  by  the  surface  regularization  algorithm.  Note  that 
while  the  data  is  being  regularized,  the  motions  are  between 
the  permitted  margin  dictated  by  the  data. 

IV.  Concluding  remarks 

In  this  paper  we  have  exploited  important  components  of 
DTED  specifications  that  have  been  widely  ignored  when 
processing  elevation  maps.  Extending  work  on  geometric 
active  surfaces  we  have  shown  how  the  maps  can  be  regu¬ 
larized  while  maintaining  the  data  in  the  intervals  allowed  by 
the  pre- specified  errors.  Additional  constraints  can  be  easily 
incorporated  to  the  framework  as  well.  The  framework  can  also 
be  applied  to  non-grid  data  properly  adapting  the  numerical 
implementation. 

We  are  currently  working  with  the  framework  described  in 
[12]  to  provide  an  alternative  to  the  approach  here  presented. 
We  are  also  investigating  the  use  of  the  results  obtained  here 
for  compression  of  elevation  maps,  [11],  and  working  on 
adding  constraints  to  the  regularization  process.  Results  in 
these  directions  will  be  reported  elsewhere. 
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Fig.  3.  Second  example  of  constrained  regularization  of  elevation  maps.  A 
set  of  original  and  regularized  iso-height  lines  are  shown,  spaced  at  16  meter 
intervals  (from  maximal  contour  elevation  of  1159  meters  to  minimal  contour 
elevation  of  919  meters). 


Fig.  4.  Quantized  vertical  adjustments  distribution.  Note  that  following  the 
imposed  constraints,  all  points  were  adjusted  less  than  10  meters  (the  cylinders 
height),  while  the  vast  majority  actually  moved  less  than  2  meters.  Other  data 
sets  might  have  a  different  distribution  for  the  points  vertical  motion  resulting 
from  the  regularization,  being  important  to  keep  the  cylinders  of  the  size 
dictated  by  the  DTED  specification,  as  done  in  this  work. 


